Surface EMG: corrections against Farina 2004, FDI layer defaults, published end-plates - #40

…bipolar

… with 1/v

The source FFT was normalised by 1/N and the inverse by N, dropping the
dz (forward) and Fs/N (inverse) factors of the continuous transforms. Their
product is dz * Fs = v, so every potential was divided by the conduction
velocity. Use dz on the forward transform and Fs on the inverse.

New test compares a single fibre against an exact, independently written
Fourier-Bessel solution for a homogeneous insulated cylinder (both routes,
two velocities, two geometries): now within 2 %, before it was off by 1/v.

The fibre source was the bare Rosenfalck slope with an ad-hoc A = 96/10 and no
fibre properties, so outputs labelled mV were not potentials. Use the published
A = 96 mV/mm^3 and the core-conductor line current density
sigma_i * pi * (d/2)^2 * d2Vm/dz2; with S/m and mm the result is in mV.

SurfaceEMG gains fiber_diameter__um (default 50) and
intracellular_conductivity__S_m (default 1.01), passed to every fibre on both
routes and included in the execution cache key. The experimental hybrid path
(use_unified=True) keeps its own D1 source and is unchanged.

Tests: exact-solution comparison now uses the physical source; amplitude scales
with d^2 and sigma_i; both routes forward the parameters.

…eld 1997)

The previous default, 4.88e-4 S/m (measured composite skin, Botelho et al.
2019, tissue data at 150 Hz), raised surface amplitudes 3.5-6x. 1 S/m is
the value used by Farina et al. (2004), taken from Roeleveld et al. (1997),
who chose it so their layered cylinder matched how biceps MUAPs fall off
with depth and around the arm. It is an effective parameter of the
cylinder model, not a measured skin conductivity. Pass 4.88e-4 explicitly
to reproduce the old default. The execution cache key already includes
the skin conductivity.

CHANGELOG: amplitude fixes, new SurfaceEMG fibre parameters, upgrade notes.

…mogeneous solution

nmj_branch_schedule, two_layer_nmj and nerve_delays__ms hold the
Konstantin et al. (2020) branch layout, end-plate positions and axon
path lengths that MotorUnitSim used inline. Draw order is unchanged, so
intramuscular EMG reproduces its previous end-plates and paths exactly.

…e grid

v_grid fixes the (k_z, k_t) grid and window; the fibre velocity enters
only the source term, whose IAP spectrum is read at k_t / v. Fibres at
different velocities then sum on one time axis. delay_fibre_signal
applies a fractional-sample delay. v_grid defaults to v, which leaves
every existing call unchanged.

…ocities

SurfaceFiberDispersion switches each effect on independently: two-layer
end-plate scatter (full or within each motor unit), axon-branch delays,
size-ordered and within-unit fibre velocity spread, and diameter-scaled
source strength. Draws use their own seeded generators; the default
(off) is bit-identical to the undispersed model. The realised spread is
reported per motor unit. The execution route and use_unified refuse an
active configuration.

…iameter line

Add fiber_cv_threshold_slope (Del Vecchio et al. 2017: motor-unit CV rises
linearly with recruitment threshold). Read fibre diameter from the
Nandedkar & Stalberg (1983) line cv = 0.95 m/s + 0.05 m/s/um * d instead of
an unsourced slope and reference diameter, and report the mean diameter.
Drop an unverified citation for the within-unit CV spread.

Konstantin et al. (2020, eq. 5-6) print a_mu = 1, b_mu = 2.5, a_sigma = 0.25,
b_sigma = 1 mm; the reference simulator code (and MyoGen's default) uses
1.5, 4, 0.5, 1.5. Offer both in nmj_branch_schedule and SurfaceFiberDispersion,
default unchanged.

Intramuscular and surface EMG now use the z spreads printed in Konstantin
et al. (2020, eq. 5-6). Within-unit end-plate SD falls from a median of
about 1.3 mm to 0.8 mm (largest units about 1.9 mm instead of 2.9 mm),
in line with simulator conventions of about 1 mm and the 2-3.5 mm upper
bounds from single-fibre and surface-array studies. The reference-code
values remain available as coefficients="reference_code".

tendon_end_band__mm places each fibre end independently, uniform over a
band about +-mean_fiber_length/2, instead of the length variation that
moves both ends together. Off by default; the dispersion realisation now
carries each fibre's ends and the surface model reads L1/L2 from them.

tendon_end_correlation__mm makes each fibre end follow a smooth Gaussian
random field over the cross-section (one per tendon, squared-exponential
covariance via random Fourier features) plus a per-fibre jitter, keeping
the total SD of the uniform band. Neighbouring fibres then end close to
each other. 0 keeps the independent band. Off by default.

…half-width

SurfaceFiberDispersion.literature_candidate() collects the values chosen
from the literature: published two-layer end-plates, 5/2 m/s nerve
branches, threshold-linear unit CV (Del Vecchio 2017), 5 % within-unit
CV and correlated tendon ends (the last two are assumptions). The
diameter-scaled source stays off. Opt-in.

IntramuscularEMG(nmj_coefficients="reference_code") restores the
end-plate spreads used up to 0.13.0; the default is the published
Konstantin et al. (2020) values.

benchmarks/semg_dispersion: paired sweeps of the SurfaceFiberDispersion
options on a limb-like cylinder and the default (FDI) muscle; arms run
before the published end-plate default pin reference_code coefficients
and reproduce the stored outputs bit for bit. benchmarks/semg_versions:
the before/after comparison of the surface-EMG fixes.

CHANGELOG: fibre dispersion, the literature preset, the published
end-plate spreads as iEMG default (nmj_coefficients restores the old
ones) and the fibre-length docstring. Coefficient wording corrected to
1.5-2x; preset docstring drops an unverified delay comparison.

gen_arm.py takes _fat_thickness__mm and _skin_conductivity__S_per_m;
existing arms still reproduce bit for bit. Eight FDI arms with skin
1.2 mm (Li 2018), fat 0.3/1.3/2.0 mm and skin conductivity 0.15 or
1 S/m, undispersed and with all effects.

The default comes from Farina et al. 2004, who took it from the
three-layer model of Roeleveld et al. 1997 (also Blok et al. 2002). It
was chosen so the model matched how biceps MUAPs fall off with depth
and around the arm; the authors note it is not known whether skin
conducts that well. Measured composite skin is ~4.6e-4 to 4.9e-4 S/m
(Gabriel 1996 via Lowery 2002 and Botelho 2019). Replaces the earlier,
unsourced 'gel bypasses the stratum corneum' rationale. Value unchanged.

Skin 1.29 -> 1.2 mm, the hand-dorsum median of healthy controls in Li et
al. 2018 (range 0.79-1.50 mm). Fat 0.3 -> 1.3 mm so that the total depth
to muscle is 2.5 mm, the low end of the FDI subcutaneous depth measured
by Nishikawa et al. 2024 (group means 2.54-3.30 mm); the fat value is
derived from that total, not measured. Default-Muscle surface amplitudes
fall by about 40 %. Skin conductivity stays 1 S/m pending further work.

On the default SurfaceEMG route any exception while simulating one motor
unit was logged and replaced by an all-zero MUAP, which passes for a
valid, silent unit. Seen with sampling_points_in_theta_domain >= 256,
where Bessel under/overflow leaves one conductor system singular. The
exception is now re-raised with the unit and array index noted; the
ExecutionPolicy route already raised.

_fiber_diameter__um sets SurfaceEMG(fiber_diameter__um); the new arms
repeat the literature-layer arms at the Hwang et al. 2013 FDI mean.

…2 rates

_rate_model selects the De Luca & Contessa (2012) exponential, with parameters
anchored to De Luca et al. (1982) Table 1 FDI rates (8.9 Hz recruitment, 25.3 Hz
peak at 40 % MVC, 41.4 Hz at 80 %). The fit is a modelling assumption: Table 1
gives per-unit peak rates of decomposed units, not the population mean.

…dent solutions

An independent Fourier-Bessel solver for an insulated cylinder with anisotropic
muscle, fat and skin (checked against the homogeneous closed form) is compared
with MyoGen at FDI layers, 0.5 and 2 mm depth, skin 1 and 0.15 S/m (1 %).
A sealed-end cable source (end-plate generation, tendon extinction) convolved
with that solver's point-source potential matches the bipolar finite-fibre
MUAP at 5 and 10 mm IED (1 %).

_muscle_radius__mm and _fiber_density__per_mm2 set the FDI cylinder;
_spike_file replaces the generated trains with NEURON spike trains for the
same threshold array (checked on load). Arms: r 5.3 mm, 267 fibres/mm²
(~40 k fibres), fat 2.3 mm, NEURON trains.

…rations

Spike files from fdi_spikes_force.py (calib/fdi-dd-deluca1982): VL-optimised drive
on the FDI pool, and FDI calibrations with recruitment ending at 70 % and 52 %
excitation, at 25 % MVC (force model) and 25 % excitation; 267 and 350 fibres/mm²,
skin 1 and 0.15 S/m. Spike-file paths point to the session scratchpad.

gen_arm.py: _delsys5 records a monopolar 3x3 grid (2.5 mm pitch, 0.25 mm
radius pins) per (z, rotation, theta) so the Delsys channels (corner minus
centre, adjacent corners) are formed afterwards; _fs__Hz sets the output rate
(10 240 Hz for the 20 Hz-2 kHz band of Botelho et al. 2019).

SurfaceEMG(iap_lambda__per_mm=1.0) sets Vm = 96 (lz)^3 e^{-lz} - 90 mV.
The default is the published 1 / mm and leaves outputs unchanged. Values
below 1 stretch the IAP and lower the MUAP's spectral content.

On the FDI defaults, 0.7 matched the median frequency and adjacent-pin RMS
of the Botelho et al. (2019) Delsys recordings at 25 % MVC (one subject);
1 / mm gives MDF ~200 vs ~165 Hz measured. Any value other than 1 is a
fitted modelling assumption; the default is not changed.

Passed through the default and execution= routes (and the disk-cache key);
use_unified=True refuses lambda != 1 because its source would ignore it.
Tests: single fibre against the exact homogeneous cylinder at lambda 0.7
and 0.5, parameter passed to every fibre on both routes, invalid values
rejected. Full suite: 389 passed, 2 skipped.

…catter

replicate.py smooths a single-fibre profile (1 mm rows along z, gen.py
z_rows_mm) with a Gaussian of SD 3 mm, which in the z-unbounded cylinder
is the expected MUAP of a point unit whose end-plates scatter by that SD.
All five replication ratios pass the frozen +-15 % limits (worst -10 %);
with SD 0 the first run is reproduced exactly. Production code unchanged.
The fall-off arms lack the scatter and still need a rerun.

…sive

gen.py places one 1 mm-row column per arc position (arcs_mm);
run_scatter_arms.sh runs skin 1 / 0.15 / 4.88e-4 S/m x territory radius
0-12.5 mm; falloff.py smooths with the 3 mm end-plate Gaussian and applies
the frozen rule. No arm meets depth 1.8 +- 0.3 and lateral <= 1.3:
territories steepen the lateral fall-off; 1 S/m is closest (lateral
1.7-2.0, depth 2.5-2.7) but still too steep; measured skin values are
further off. Production code unchanged.

Only ~100 S/m reproduces the measured lateral and depth spread; 1 S/m is
the least-bad physical-range value and measured skin values are 3-6x
further off.

MyoGen cited Akhmadeev, Yu, Le Carpentier, Aoustin & Farina (2020, IEEE
TBME 67:2005) as "Konstantin et al."; Konstantin is the first author's
given name. Docstrings, comments, the recruitment-threshold example and
CHANGELOG corrected. The konstantin threshold mode keeps its name.

The Roeleveld 1997 benchmark was inconclusive and is not evidence for
1 S/m. State the actual reason (a measured 0.15 S/m puts FDI plateau RMS
3-4x above published ranges) in the CHANGELOG and Muscle docstring.
Move the upgrade notes to a Breaking changes section and record the
use_unified bipolar and angular-sampling fixes.

…ace-EMG fixes

arms.json pointed at spike files in a local temporary folder. The five
files these arms use (120 units each, ~2 MB) now live under
benchmarks/semg_dispersion/spikes/ and _spike_file paths are relative to
the benchmark directory.

The calibration argument for 1 S/m compared against a 'measured 0.15 S/m'
without a source, next to the 4.88e-4 S/m composite skin also called measured.
0.15 S/m is the IT'IS tissue database's low-frequency skin value (0.148 +/-
0.042 S/m, 7 studies, 0.057-0.25).